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Abstract 


Around a million quasars have been catalogued in the Universe by probing 
deeper and using new methods for discovery. However, the hardest ones to find 
seem to be the rarest and brightest specimen. In this work, we study the properties 
of the most luminous of all quasars found so far. It has been overlooked until 
recently, which demonstrates that modern all-sky surveys have much to reveal. The 
black hole in this quasar accretes around one solar mass per day onto an existing 
mass of ~17 billion solar masses. In this process its accretion disc alone releases 
a radiative energy of 2 x 10%! Watts. If the quasar is not strongly gravitationally 
lensed, then its broad line region (BLR) is expected to have the largest physical 
and angular diameter occurring in the Universe, and will allow the Very Large 
Telescope Interferometer to image its rotation and measure its black hole mass 
directly. This will be an important test for BLR size-luminosity relations, whose 
extrapolation has underpinned common black-hole mass estimates at high redshift. 


In 1963, Maarten Schmidt identified the first quasar! known as 3C 273. It appeared 
as a remarkably bright star of 12th magnitude, while its redshift suggested that it was 
among the most distant objects known in the Universe at the time. The two facts 
together implied an implausibly huge output of light, and ever since then, newly found 
quasars have impressed with their immense energy release from a small region of space. 
This could only be explained as gravitational energy being converted into heat and light 
within a highly viscous accretion disc around a supermassive black hole2Bl (SMBH). 
Quasars are signposts of fast growth in SMBHs on public display and allow the study 
of these growth processes. 

Finding large samples of quasars then provides population and growth statistics to 
explain the origin of SMBHs in the Universe Generally, the most luminous quasars 
contain the fastest-growing SMBHs, although the relation between mass accretion rate 
and luminosity is affected by the mass and spin of the black hole as well as the structure 
and viewing angle of the accretion disc and disc winds BEA 

Today, around a million quasars are known, although some specimens stand out 
from the crowd: in 2015, the ultra-luminous quasar J0100+2802 at redshift z = 6.3 
was identified“© with a SMBH of 10 billion solar masse In 2018, an even more 
luminous object was found, J2157—3602 at z = 4.7, with a SMBH of 24 billion 
solar masses”! Although their luminosity implies rapid growth, their existence is hard 
to explain: when black holes start from the remnant of a stellar collapse and grow 
episodically within the Eddington limit, they are not expected to reach the evident 
masses in the time from the Big Bang to the epoch of their observation, which has 
triggered a search for alternative scenarios 2122] 

While exceptionally rare, the most extremely luminous quasars are interesting for 
several reasons beyond their intrinsic nature as discussed later. In this paper, we present 
the properties of the recently discovered? quasar SMSS J052915.80—435152.0, here- 
after J0529—4351, which is a 16th magnitude object at redshift z = 3.962, and reveal 
it to be the most luminous quasar currently known in the Universe (see Fig.[I). 


Results 


When quasars appear extremely bright, it may be suspected that their observed bright- 
ness is magnified by gravitational lensing from a massive galaxy on the line of sight. 
Strong lensing causes multiple separate images of a quasar in the sky2423] Two other 
quasars with a redshift and apparent brightness similar to JO529-4351 are known to be 
strongly magnified by lensing, the double-image APM 08279+5255 at z = 3.91, with 
a separation of 0.5 arcsed24] and the quadruply imaged B 1422+231 at z = 3.62, with 
separations up to 1.5 arcsec“. Estimated magnification factors for these two objects 
range from 40 to 100, which implies that these quasars are not intrinsically extreme, 
but are members of the bulk population2628) 

JO529-4351, in contrast, shows no sign of strong lensing; data from the European 
Space Agency (ESA) Gaia satellite suggest it to be a point source, in terms of object 
morphology and astrometric excess noise, which has been used to find dual quasars 
or lensed quasars that appear unresolved to Gaid? (see Figure [1]and online methods 
for more detail). We also find no strong foreground absorber system, which our high 


signal-to-noise spectrum probes in MgII at z > 1.15; the strongest system, at z = 
2.118, has an equivalent width of EW(2796)=0.8A; this suggests an impact parameter 
of 20 kpc (or 2.4 arcsec) from the line of sight to the quasar, while a plausible image 
separation in a lensed scenario is ~ 0.1 arcsec. We can also estimate the probability 
of lensing a source at z ~ 4 using common models for the galaxy distribution as 
isothermal masses3132) obtaining p ~ 1.3 x 107°. An image separation of 0.2 arcsec 
or less reduces it further to p ~ 2 x 1074. These estimates change by only a factor 
of a few when changing the input galaxy velocity dispersion function or the bright 
end of the quasar luminosity function. Even if a very steep intrinsic quasar luminosity 
function at the bright end enhanced the magnification bias, the probability that this 
source is strongly lensed will be less than 1%. We thus take the strength of the quasar 
emission at face value, although final confirmation from a high-resolution space-based 
or adaptive-optics image would still be desirable. 

We also investigate the recent history in the brightness of JO529-4351 to see whether 
it may have been previously overlooked due to extreme variability. The 0.5-metre 
NASA (National Aeronautics and Space Administration) Asteroid Terrestrial-impact 
Last Alert System?! (ATLAS) telescope provides a light curve since mid-2017 (see 
Fig. B). Brightness variations of ~ 15% are found over the last six years, which are not 
unexpected for luminous quasars. Earlier records, from photographic plates observed 
in 1980 and 1998 (measured by the SuperCOSMOS Sky Surveys), found the R-band 
brightness to be consistent with recent observations by the SkyMapper Southern Sur- 
vey! (SMSS). The somewhat longer light curve from the Wide-field Infrared Survey 
Explorer (WISE) shows also only modest variability. 

The quasar is undetected in the Rapid ASKAP Continuum Survey24l (RACS), and 
hence has a flux of less than 14Jy (4c) in the broad 887 Mhz band. Using a common 
definition of radio loudness we find the object to be safely in the radio-quiet regime 
(R < 1). Thus, we have no reason to suspect that its luminosity is affected by jet 
emission, let alone relativistically boosted. 

For further analysis, we use a spectrum of optical and near-infrared light (see Fig.|3) 
from the 8.2m Very Large Telescope (VLT) at the European Southern Observatory 
(ESO). We split the spectrum with a publicly available spectrum-fitting code! into 
an accretion disc continuum and emission-line contributions. From the disc contin- 
uum we quantify the disc luminosity of the quasar and get a first proxy for the ac- 
cretion rate. We find monochromatic luminosities of log(L135/erg s~') = 47.93 and 
log(L3o0/erg s71) = 47.76 at restframe wavelengths of A = 135 and 300 nm, re- 
spectively. An approximate estimate of the full radiative output from the accretion 
disc using standard bolometric correction? and a 0.75 anisotropy factor“ yields 
log(Lpoi/erg s~+) = 48.37. Using a standard radiative efficiency valud of 0.1, this 
translates to an accretion rate of ~ 413 solar masses per year. This result makes JO529— 
4351 the most luminous quasar and by inference, the fastest-growing black hole in the 
Universe known to date (in terms of mass growth per unit time). 

A more refined analysis of the bolometric luminosity and accretion rate requires 
modelling the possible spectra of the disc continuum over a grid of black-hole masses 
and spins and simulating their observation for different accretion rates and viewing 
angles 2/45) Using slim-disc models as appropriate for SMBHs with high accretion 
ratek and a publicly available code (see online methods), we find a best-fitting 


solution to our spectrum being an Eddington-accreting disc around an SMBH with a 
moderate spin of a ~ 0.4 viewed at an intermediate angle of 2 ~ 45° (see Fig. (4). 
The intrinsically emitted bolometric luminosity is log(Lpo1/erg s71) = 48.16, which 
is 62% of the value derived above; the lower value mostly results from generic bolo- 
metric corrections being overestimates for very massive black holes with their colder 
discs and suppressed UV emission. The best-fit slim disc model yields an accretion 
rate of 370 Mo per year and a radiative efficiency of ~ 0.09. Given the lack of higher- 
frequency data in the SED, the integrated model SED may still underpredict the bolo- 
metric flux, while the standard bolometric corrections must be overcorrecting for black 
holes of over a billion solar masses. The best estimate is expected to lie between the 
two, and we thus choose the average of log(Lyo1/erg s71) = 48.27. 

When modelling SEDs for a fixed mass accretion rate, different black-hole spins 
lead to different bolometric corrections and radiative efficiencies, while the monochro- 
matic UV-optical luminosities are only modestly affected. The main uncertainty arises 
from the unknown viewing angle, which affects both the apparent monochromatic lu- 
minosity and the derived accretion rate. Given the broad confidence intervals for spin 
and viewing angle in the model fit, we consider the full plausible angle range from pole- 
on to ¿ = 60° as a 95% error margin; this corresponds to a 95% range in luminosity of 
+0.12 dex and in accretion rate of 280 to 490 solar masses per year. 


Black-hole mass 


We estimate the mass of the black hole powering this quasar using two fundamentally 
different methods: 

(1) Assuming that the continuum emission is affected at its blue end by the inner 
truncation of the accretion disc due to the innermost stable circular orbit around the 
black hole! we can infer a combined estimate for mass and spin of the black hole. 
More massive black holes impose larger truncation radii and move the peak of the con- 
tinuum emission to cooler temperatures and longer wavelengths“ With this method 
the mass of the black hole in 3C 273 was found to be between 200 and 500 million 
solar masses 35 years before a value of ~300 million solar masses was measured 
by interferometric observations! of the broad-emission line region in 3C 273. Since 
then, the a o has not only proven useful for estimating masses 
using thin disc models“! but has also been applied/4l with slim discs models“ that 
are expected to be a more realistic description of the near-Eddington accretion discs 
of fast-growing SMBH. From the continuum shape of JO529-4351, we find a mass of 
log M/Mọo = 10.2871) (or ~ 19 billion solar masses). 

(2) Assuming that the broad emission-line region represents virialised gas moving 
at the velocity of Keplerian orbits around the SMBH, we can use the width of emission 
lines and the continuum luminosity to infer the mass of the black hole. This method 
is known as the virial single-epoch method=253l and has also been used to estimate the 
mass of the SMBH in JO100+2802 at z = 6.3. Its application at the high-luminosity 
end of the quasar population is mainly limited by a lack of independent calibrations for 
the SMBH mass and relies on extrapolations from lower-luminosity quasars calibrated 
with reverberation mapping=4) 


Our spectra offer two emission lines for this method, the triple-ionised carbon line 
C1vVA1549 and the singly-ionised magnesium line Mg 11\2799. The CIVA1549 line 
appears asymmetric and blueshifted relative to the Mg 11\2799 line (Av = —3120 + 
80 km/s) and a line full-width at half maximum (FWHM) of 7245 + 175 km/s. The 
Mg 11A2799 line turns out to be hard to calibrate at this redshift because of atmospheric 
absorption, but we do measure a FWHM of 4, 395+ 435 km/s. Using calibrations com- 
monly applied to high-luminosity quasars, the line properties and luminosity translate 
into SMBH mass estimates from log M/Mo = 10.03 + 0.06 to 10.45 + 0.02 for the 
CIV line=55) and from log M/Mo = 10.030.09 to 10.36+0.09 for the Mgl lineS954 
(errors are standard deviations due to propagated observational uncertainties, and the 
scatter among the values is in line with large systematic calibration uncertainties). 

The line-based estimates are consistent with each other and with the continuum- 
based estimate. Given that the systematic uncertainties in these methods are larger 
than the statistical error propagation and may be as high as 0.4 dex, we combine the 
mean mass estimates of the two lines and that from the SED without weighting into 
a final result for the mass of log M/Mọo = 10.24 + 0.02 (simple mean and standard 
deviation), or ~ 17 billion solar masses. This also implies that the SMBH is accreting 
near the Eddington limit (Eddington ratio of ~ 0.9). 


Discussion 


In terms of luminosity and likely growth rate, JO529-4351 is the most extreme quasar 
known. The accretion of JO529-4351 is near the Eddington limit, which is common 
among quasars of the highest luminosity B859 The growth rates are mostly uncertain 
due to the unknown viewing angle. Assuming persistent accretion at its current Edding- 
ton ratio, the mass doubling time is ~ 30 Myr. However, with 19 billion solar masses, 
the black hole in J0529-4351 at z = 3.962 is not the largest SMBH found in the most 
luminous quasars. It has over 50% more mass than the black hole in JO100+2802 at 
z = 6.3 but one third less mass than that in J2157-3602 at z = 4.7. The age of the 
Universe at these three redshifts is 859 Myr, 1,244 Myr and 1,530 Myr, assuming a 
flat Universe with a concordance cosmology (a cosmological constant of Qa = 0.3). 
Given that JO529-4351 is observed at a later epoch in the Universe than J2157-4351 
and JO100+2802, it is less of a challenge to models of early SMBH growth. 

It may be tempting to speculate on the spin of extreme SMBH: while it has been 
argued that growing black holes should spin up over time, this would also increase 
radiative feedback and slow down accretion, making it harder to grow the most massive 
SMBH within the age of the Universe; instead, one way for growing black holes from 
stellar seeds to the greatest masses we measure is “chaotic accretion” with randomly 
changing orientation that keeps the black-hole spins and radiative feedback low 

A long-standing question has been what mechanism fuels the high accretion rate, 
which also must have persisted for some time already, though not necessarily in the 
immediate past. Mass and kinematics of gas in the host galaxy of J0529-4351 could be 
observed with the Atacama Large Millimetre Array (ALMA). ALMA has already re- 
vealed the largest spiral galaxy in the early Universe in the quasar BRI 1335-0417 at 
z = 4.4, which is ~ 10x less luminous than JO529-4351. In contrast, ALMA observa- 


tions have revealed a merger signature in the host galaxy of a quasar at z = 7.54 Fur- 
ther observations of extremely luminous quasars are progressing? If extreme quasars 
were caused by unusual host galaxy gas flows, ALMA should see this; if nothing un- 
usual was found in the host gas, this would sharpen the well known puzzle of how to 
sustain high accretion rates for long enough to form such extreme SMBHs. 


Have we overlooked still more extreme quasars? 


Finding rare and exceptionally bright quasars such as J0529-4351 does not require large 
telescopes, but is instead a needle-in-the-haystack problem that needs precise data with 
discriminative power across large areas of sky. Such objects are often hiding in plain 
sight and are mistaken for stars. E.g., the quasar J1144-4308 at redshift z = 0.83 was 
recognised as one of the brightest quasars in the night sky (at 14 mag) only in 2022, 
despite having been imaged on photographic plates since the 19th century®31 

In the 1960s finding the first quasars was driven by radio detections of objects that 
otherwise appeared like regular nearby stars in our Milky Way Galaxy. i ea 
surveys employed colour selection, template-fitting and Bayesian methods“2! But 
the main barrier for obtaining complete samples has always been contamination of 
quasar candidate selections with stars from our Milky Way, which appear similar if 
no discriminating information is available and vastly outnumber true quasars in the 
extreme regime of the brightest quasars. With a million quasars known by now and 
the wealth of data from modern all-sky surveys, machine-learning approaches are now 
most popula E5] These tend to get the classification right for the majority of objects, 
while training samples imply that they perform less well for rare extreme cases. In the 
case of extremely luminous quasars, the obvious bias of a training sample is simply 
that they do not seem to exist until they are found. 

A low-resolution spectrum of J0529-4351 revealing its quasar nature and redshift 
has been part of the public all-sky Gaia DR3 data set published on 13 June 2022. The 
machine-learning classification of this data set by the Gaia Discrete Source Classifier 
(DSC) has assigned J0529-4351 a 99.98% probability to be a Milky Way sta, al- 
though a human astronomer eyeballing the Gaia spectrum would recognise the quasar 
and redshift at first sight. 

These days, combining data from Gaia and WISE makes an all-sky search for bright 
quasars straight-forward. An absence of parallax and proper motions removes most 
of the bright Galactic stars from the search, and the mid-infrared photometry from 
WISE is then sufficient to discriminate the disk continuum and dust emission of quasars 
from the Rayleigh-Jeans tail of stellar photospheres® The All-sky Bright Complete 
Quasar Survey=!(AIIBRICQS), which reported JO529-4351, selects its candidates with 
a simple heuristic WISE colour cut of W1—W2> 0.3 and an absence of > 4c evidence 
for stellar motion. While it aims for completeness and its recall of known quasars is 
nearly perfect, its efficiency is still high, given that 97% of its candidates turn out to 
be quasars) Given that the AIIBRICQS follow-up is complete at Rp < 16 in both 
hemispheres, we doubt that a quasar of higher UV-optical luminosity will be found in 
the future, unless it is hiding behind the Galactic Plane: AIBRICQS did not search at 
Galactic latitude of |b| < 10°, where source crowding and higher dust extinction makes 


the discovery of quasars still very challenging) 


Outlook: future observations with bright quasars 


Extremely luminous quasars enable further observations for specific challenging quests: 
after decades of instrumental development to increase the spatial resolution of tele- 
scopic observations, the VLTI/GRAVITY Collaboration= made headlines with their 
spatio-kinematic mapping of the broad-line region (BLR) in the iconic nearby quasar 
3C 273 (redshift z = 0.157). This observation revealed the orbital rotation of the 
disk-like BLR in 3C 273 in a spatially resolved pattern and thus provided a direct mea- 
surement of the black hole mass from the BLR orbits. The measurement was possible 
as 3C 273 was bright enough and displayed its Paschen-a emission line in GRAVITY’s 
K-band window. 

While the redshift of JO529-4351 imposed challenges for the work presented here, 
as nearly all strong UV emission lines are in places heavily affected by atmospheric 
absorption, it does place its H£ line at 2.44. in the K-band window of the GRAVITY 
instrument. We predict the size of the HG-emitting broad-line region in JOS29-4351 by 
extrapolating the radius-luminosity relation®l and find a radius of 2.2 pe. This implies 
an angular diameter of 0.64 milli-arcseconds, which is an order of magnitude larger 
than the BLR in 3C 273 and thus the largest-appearing quasar BLR in the Universe. We 
expect the soon-to-be-upgraded VLTI/GRAVITY+ to obtain a superbly well-resolved 
picture of the BLR rotation in JO529-4351 and thus a much more reliable measurement 
of its black hole mass. Crucially, black-hole masses at the high end are estimated from 
relations that have been extrapolated by orders of magnitude such that the whole scale 
for objects like this is at risk. Getting a direct mass measurement for a black hole 
with likely 50x the mass of the black hole in 3C 273 would be extremely valuable 
for constraining the relations commonly used to estimate the masses of black holes in 
the early Universe. The true mass scale of the earliest SMBH would also impact the 
question of how hard exactly it is to form them. 

Future plans include watching the Universe expand with repeat observations of 
quasars lasting for a decade: the expansion shifts the redshift of individual gas clouds 
observed as Lyman-a forest absorption lines in quasar spectra 281] (Sandage test). The 
signal of such observing campaigns depends on the availability of a sufficient number 
of very bright quasars and will still require advanced facilities such as the forthcom- 
ing ESO Extremely Large Telescope (ELT). JO529-4351 will quite obviously be an 
important part of this long-term endeavour. 


Methods 


Evidence for point-source geometry 


While Gaia provides no images as such, we use the source characterisation from 
its Data Release 3 to provide evidence for a point-source geometry in J0529-4351. 
First, there is no noticeable flux outside an aperture with 0.175 arcsec radius around 
the object centroid in the high-resolution imaging®2) as determined by comparing the 
fluxes between Bp/Rp and the G-band aperture) The closest neighbouring object 
is 3.2 arcsec away; it is not detected by Gaia, but by the DECam component of the 
DESI Legacy Imaging Survey $? with g = 23.581 + 0.039, r = 22.319 + 0.015, 
i = 20.686 + 0.007, and z = 20.099 + 0.007 in DR10, which is a good match to the 
colours of an M4 star in the Milky Way. 

Second, astrometric excess noise is absent in J0529-4351, although 319 good ob- 
servations were available to determine astrometric solutions. All Gaia parameters that 
may hint at unresolved multiple sources, including frac_multipeak (0), 
astrometric_excess_noise (0), astrometric_sigma5d_max (0.069), 
ruwe (1.013), and ipd_gof_harmonic_amplitude (0.0167) are in line with the 
average values (0.132, 0.042, 0.074, 1.015, 0.045) for bright quasars with R, < 16 at 
redshift z > 2 and far below any selection cuts used to search for Gaia-unresolved 
binary stars®4l and dual or lensed quasars Figure 1 shows the astrometric excess 
noise (AEN) of quasars at R, < 18 including all dual and lensed quasars2285] with im- 
age separations below 0.5 arcsec from the list athhttps://research.ast.cam.| 
that are unresolved by Gaia; even at the smallest known 
image separation of 0.18 they show an (AEN)> 1.0 mas. 

For the lensing calculation we followed a procedure used in the analysis in a re- 
cently discovered high-redshift lensed quasar? while adjusting the source redshift to 
the z ~ 4 of JO529-4351. 


Spectroscopic observations and data processing 


The eee spectrum of JO529-4351 was obtained with the X-Shooter instru- 
menton ESO’s Very Large Telescope (UT3) on 2023, Jan 14. The airmass at the time 
of observation was 1.06 and the estimated seeing was about 1.5 arcsec. In the UVB, 
VIS, and NIR arms, the slit widths and total exposure times were: 1.0 arcsec and 500 s; 
0.9 arcsec and 400 s; and 0.9 arcsec and 800 s, respectively. 

The data were reduced using Pypelt, a Python-based spectroscopic data reduction 
pipeline®288| Because it was not possible to a priori know the underlying true shape 
of the C 1V\1549 emission line beneath the atmospheric A-band absorption, we could 
not reliably fit the telluric model to the JO529-4351 VIS spectrum alone. To obtain 
an improved correction, we took the telluric model generated from a standard star 
observed on the same night and modified it by adopting a power-law intrinsic stellar 
spectrum across a broad span of wavelengths unaffected by strong stellar atmospheric 
lines. Applying this modified model to JO529-4351 yielded a final spectrum with a 
smooth CIVA1549 profile containing only narrow absorption lines and no spurious 
emission from overcorrected telluric absorption. 


The observed spectrum was corrected for slit losses by calibration to external pho- 
tometric data, in particular, the quasar’s VHS DR6 J-band magnitude of 14.812 + 
0.003 mag (Vega). With the VHS calibration anchoring the NIR arm, the VIS and 
UVB data were sequentially matched to the spectroscopic flux levels in the wavelength 
regions of overlap between the arms. This J-band calibration is consistent (within the 
photometric errors) with that of the zgpss X-Shooter acquisition image, as calibrated 
by synthetic photometry for a neighbouring, non-variable star with Gaia low-resolution 
spectroscopy, and the flux scale measured in the Gaia spectrum of J0529—4351 itself. 


Spectral decomposition and emission-line fitting 


The reduced and telluric-corrected X-Shooter spectrum is then transformed into the 
rest-frame using a redshift of z = 3.962. We model the broad emission-line profile 
of both the CIVA1549 line and Mg 11\2799 line using the code PyQSpecF it Bl a 
python-based spectral modelling package that is designed for rest-frame UV and op- 
tical quasar spectra. We model each line individually, selecting independent windows 
in wavelength to constrain the quasar continuum and the emission-line flux. We mask 
the narrow absorption features present in our spectrum by applying a boxcar sigma- 
clipping routine with a width of 50 pixels and a 3ø threshold. 

For CIVA1549, we constrain the continuum on either side of the feature with a 
power-law model over line-free wavelength windows of 1445A-1455A and 1973A- 
1983A. The emission-line model is composed of a maximum of three broad Gaussian 
components with a minimum full-width at half-maximum (FWHM) of 1000 km s~!. 
We fit the CIVA1549 profile between 1480A-1528A and 1537A-1565A, avoiding the 
narrow absorption feature between 1528A-1537A. We also measure the monochro- 
matic luminosity at 1450A using the power-law model of the continuum. We estimate 
the errors in the spectral fitting by creating 50 realisations of the spectrum after ran- 
domly redistributing the flux data according to each point’s Gaussian uncertainties. The 
standard deviation in each of the derived properties of the C 1VA1549 emission-line is 
adopted as the statistical error. 

For Mg 11A2799, the continuum model is composed of a power-law and a template 
of the UV flux contribution from broad Fe II emission, which is constrained across 
the wavelength ranges 1973A-1983A, 2060A-2340A, 2600A-2740A, 2840A-3100A. 
The emission-line model for Mg 11å2799 is composed of up to one narrow and three 
broad Gaussian components, where the narrow and broad distinction is set at a FWHM 
of 1000 km s~!. We fit the Mg 11,2799 line between 2700A-2870A and the monochro- 
matic luminosity as 3000A is also measured from the continuum power-law model. We 
derive the uncertainty in each of the line properties by adopting four different templates 
of the UV Fe 11 emission®2&2| As with C1vA1549, we also resample the spectrum 50 
times and measure the standard deviation in the Mg I1A2799 line properties, summing 
in quadrature with the uncertainty derived from Fe II templates. More information 
about the PyOSpecFit modelling of CIVA1549 and Mg IIA2799 lines are detailed 
in a study of the now second-most luminous quasar, SMSS J2157-3602441 and also 
in a study of line properties in the luminous quasar sample of the European Southern 
Observatory Large Programme XQ-10052, which is at a comparable redshift range to 
JO529-4351. 


Luminosity calculations 


For calculating monochromatic continuum luminosities, we corrected the spectra for 
extinction by dust in the Milky Way. We used the estimate of E(B — V) = 0.041 + 
0.0027 from the Schlegel maps”) with the correction factor™ of 0.86 and the extinction 
law by Fitzpatrick®) which yields restframe absorption values of A1450 = 0.074 mag 
and A3000 = 0.025 mag. We then apply bolometric corrections for mean quasar 
SEDs22! of kı350 = 3.81 and k3999 = 5.15 and apply an anisotropy correction factor 
of 0.75 assuming mean orientation“2). The accretion disc continuum fitting described 
below determines the bolometric luminosity directly from the integrated disc model 
SEDs and implies a spin- and inclination-dependent anisotropy correction. We then 
take the average of the luminosity measurements. 


Spectral decomposition and continuum fitting 


We use a publicly available code“8l to model the shape of the accretion disc continuum 
to spectral energy distribution models predicted for slim discs with s Limbh#d which 
is a grid of synthetic spectra from ray-traced numerical solutions of slim accretion disc 
equations. The free parameters are mass and spin of the black hole as well as accre- 
tion rate and disc inclination. A Markov chain Monte Carlo (MCMC) method is used 
to map out likelihood contours across the 4D parameter space. We collect the Milky 
Way extinction-corrected spectrum and extend the IR coverage with WISE photometry, 
specifically the CatWISE202024 W7 and W2 photometric passbands. We then create 
synthetic data points to represent the accretion disc flux within selected line-free win- 
dows in the observed spectrum. For the WISE wavelengths, we adopt the spectrum 
of the Selsing high-luminosity quasar template) scaled independently to the flux of 
W1 and W2 bandpasses, and create synthetic line-free data points in the same fashion. 
Furthermore, we use the X-Shooter UVB arm to estimate the continuum flux through 
the Lya forest and create one more synthetic data point to constrain the continuum. We 
then use MCMC to infer the Bayesian posterior probability distributions of the intrinsic 
black hole properties using the set of synthetic data points to represent the flux of the 
accretion disc continuum. 

Further detailed information on the continuum fitting method is presented in an in- 
depth study of the quasar SMSS J2157-360241 However, unlike this previous study, 
the static and smooth thermal accretion disc models are unable to fully reproduce the 
hardening of the continuum in the VIS arm or the flux passing through the Lya forest. 
It has been shown previously that quasars with high C Iv\1549 blueshift, which are 
indicative of strong outflows and winds, tend towards a bluer UV continuum to the 
blue side of Arest 200 nm and a slightly redder continuum on the red side] thus 
producing a spectral break that is missing in the thermal disc models. Therefore, our 
slimbh synthetic spectra utilises the disc atmosphere model, BHSPEC2M to help 
reproduce the emerging Compton-hardened radiation. 

We also a to fit the data with alternative models. E.g., we adopt kerrbb 
thin disc models! despite the fact that the slim disc model reproduces the thin disc 
SED at low Eddington ratios, and find indeed similar results for the mass of the black 
hole and the luminosity of the disc, although the colour of the UV continuum is not 
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properly reproduced in the absence of Comptonisation considerations. We also exam- 
ine the option of JO529-4351 being lensed by demagnifying the spectrum and find that 
it does not improve the quality of fit. 


Host galaxy dust extinction 


We assume no dust extinction in the quasar host galaxy given that typical host redden- 
ing levels are found to be consistent with E(B— V) ~ 0.0 among luminous quasar [A 
while less than 1% of quasars seem to have E(B—V) > 0.1. Also, an object appearing 
as the brightest object in the Universe has a low probability of being extinguished by 
notable levels of dust. If dust were present, it would make the continuum redder and 
the measured luminosity lower. This will lead to the black hole mass being overesti- 
mated by the continuum-fitting method and underestimated by the virial single-epoch 


method. Such a discrepancy is not observed here. 


Light-curve construction 


The light curve with photometry from the NASA (National Aeronautics and Space 
Administration) Asteroid Terrestrial-impact Last Alert System] (ATLAS) 0.5m tele- 
scope was obtained from the ATLAS website. For any object in the ATLAS observing 
area, up to four observations per night are available, depending on weather. The orange 
passband is observed in all clear nights, while the cyan passband is only used during 
the half period around New Moon. For slowly varying objects, noise can be suppressed 
by combining observations from longer periods; we determine median measurements 
to reduce the influence of outlier measurements, using per-week intervals in the orange 
passband and per-moon period intervals in the cyan passband. Error bars on the medi- 
ans express the level of variability within the interval by showing inter-quartile ranges 
of the values. 


Data availability 


Data from Data Release 3 (DR3) of the European Space Agency’s Gaia mission are 
publicly available (https://gea.esac.esa.int/archive/). NASA ATLAS data are available 
from The SkyMapper Southern Survey data 
are available from (doi: 10.25914/5f14eded2d116). The 
reduced spectrum is available from the authors on reasonable request. 


Code availability 


The spectral fitting code and the quasar continuum fitting code were written by SL in 
python and are i o available on ; ithub at|https://github.com/samlaihei/PyQSpecFit=! 


and https:// github.com/samlaihei/BADFi FS) 
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Table 1: Luminosity, black hole mass and accretion rate for the quasar JO52915.80- 
435152.0: given are best-fit estimates and 68% confidence intervals. Systematic uncer- 
tainties for the mass estimates are proposed to be as high as 0.4 dex. 


Luminosity log L/(erg s71) 
AL 135 47.93 

AL 300 47.76 

Lho (from 135 nm) 48.39 

Lho (from 300 nm) 48.35 

Lro (slim disc model) 48.16 

Lro (best estimate) 48.27 + 0.06 
Black hole mass log M/Mo 
CIvA1549 line (average) 10.24+0.15 
Mg I1\2799 line (average) 10.20 + 0.08 
Continuum fit (slim disc 10.28+2-17 

All methods (best estimate) 10.24 + 0.02 
Mass accretion rate log M/(Mo yr—') 


for plausible inclination range 2.57 + 0.06 
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Figure 1: Left: The quasar JO529-4351 and other quasars that are remarkably luminous 
for their cosmic epoch in relation to the population of known quasars from SDSS and 
from XQ-100, which has the highest-quality spectra2! for luminous quasars at z ~ 4. 
The luminosity axis shows monochromatic luminosity with standard bolometric cor- 
rections applied. Bottom left: Astrometric excess noise (AEN) of J0O529-4351 com- 
pared to the population of known bright quasars at z > 1; the known Gaia-unresolved 
lensed or dual quasars are all at AEN> 1 mas. Zero AEN was replaced with a value of 
0.01. Bottom right: Image from the Dark Energy Camera Legacy Survey DR10 with 
size 20 x 20 arcsec”; North is up and East is to the left; a neighbouring M star appears 
in red. 
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Figure 2: Light curves of JO529-4351: data from NASA ATLAS are combined into 
weekly (orange passband) or moon-period (cyan passband) median averages. Data 
from WISE are median averages of multiple visit completed within less than two 
weeks. Error bars show interquartile ranges within the averaging period. 
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Figure 3: Spectrum of the JO529-4351 from the ESO 8.2m Very Large Telescope. Top: 
observed-frame wavelength range of 0.42 — 2.27um; a power law with a, = 0.5 is 
shown as a dashed line for comparison. Centre: decomposition, emission-line fits, and 
noise level (grey line, very low in the C IVA1549 part): derived black hole masses are 
log M/Mọo = 10.24 for CIVXA1549 (left) and 10.20 for Mg 11\2799 (right). Bottom: 
earlier optical-only spectra from AIIBRICQS and from Gaia DR3. 
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Figure 4: Results of the continuum fitting for JO529-4351 using MCMC and slim 
disc models. The highest likelihood model is a disc viewed at ~45° angle around 
an intermediate-spin SMBH with ~19 billion solar masses and an accretion rate of 
370 Mọ per year; the constraints around viewing angle and spin are, however, weak. 
Spin is denoted by a and Laisc is in units of Eddington luminosity. 
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